CTRW Pathways to the Fractional Diffusion Equation. 
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Abstract 

The foundations of the fractional diffusion equation are in- 
vestigated based on coupled and decoupled continuous time 
random walks (CTRW). For this aim we find an exact solu- 
tion of the decoupled CTRW, in terms of an infinite sum of 
stable probability densities. This exact solution is then used 
to understand the meaning and domain of validity of the frac- 
tional diffusion equation. An interesting behavior is discussed 
for coupled memories (i.e.. Levy walks). The moments of the 
random walk exhibit strong anomalous diffusion, indicating 
(in a naive way) the breakdown of simple scaling behavior 
and hence of the fractional approximation. Still the Green 
function P{x,t) is described well by the fractional diffusion 
equation, in the long time limit. 



I. INTRODUCTION 

Fractional calculus is an old field of mathematical anal- 
ysis which deals with integrals and derivatives of arbi- 
trary order [0-^. Fractional diffusion equations were 
introduced to describe anomalous non Gaussian trans- 
port systemsJ0-O] (see |Q for review). The stochastic 
foundation [|l4|-|l7| of these equations is the continuous 
time random walk (CTRW) introduced by MontroU and 
The relation between CTRW and the frac- 
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tional equations is the reason for a renewed interest in 
the properties of CTRWs. 

In this work we investigate the limitations and the do- 
main of validity of the fractional diffusion equation based 
on coupled and decoupled CTRWs. Some limitations 
on the fractional framework were partially addressed in 
[ p^ , p7|Jl9| , for the sub diffusive case in dimension d. Here 
we consider the one dimensional fractional diffusion equa- 
tion il 



9" 



' P,.(.,t) + 5#i^, (1) 



r(i - a) 



where is the fractional Riemann-Lioville (time) 

derivative of order a and is the Riesz space frac- 

tional derivative of order /i. These fractional deriva- 
tives are integro-differential operators, whose definition 
is given in j^,^,|l3| . The last term on the right hand side 
of Eq. (|l|) is the source term which depends on initial con- 
ditions. We consider free boundary conditions and initial 



conditions concentrated on the origin Pfr{x,0) = S{x), 
then the Fourier-Laplace (k-u) transform of the Green 
function is 



Pfr{k,u) = 



\k\f^/2' 



(2) 



This equation is for our purposes the definition of the 
fractional equation (mj . The inversion of Eq. (H) yields 



(3) 



where K{z) is a scaling function whose properties are 
given in (|J2^ . The probability interpretation of Pfr [x, t) 
is restricted to {0 < ^ < 2} n {0 < a < 1} and 
{1 < a < ^ < 2} When ^ = 2 and a = 1 the 

fractional equation reduces to the ordinary Gaussian dif- 
fusion equation, if /i < 2 and a = 1 it describes Levy 
flights. When fi — 2 the equation describes sub or en- 
hanced diffusions, (x^) oc t", according to a < 1 or a > 1, 
respectively. 

Eq. (^ has a long history: for certain values of a and fi 
it was derived from the CTRW model §Jl|,|2l|-|3 , using 
the long wave length A: — > small u — >■ approximation 
(see details below). This approximation was used many 
times to investigate the long time behavior of the CTRW. 
It is based on the simplifying assumption that the scal- 
ing behavior Eq. (||) holds |2^ . In 1 25 2^ a rigorous 
approach, based on limit theorems, was used to classify 
the asymptotic behaviors of different types of CTRWs. 
The work in is in agreement with Eq. (||) and the 

older work in this field (i.e., again for certain a and fi 
and see details below). 

Here exact solution of the decoupled CTRW in {x, t) 
space is found in terms of an infinite sum of stable proba- 
bility densities. This exact solution is used to investigate 
the meaning and limitations of the fractional diffusion 
equation. For example: we show that certain solutions of 
the fractional diffusion equation diverge on the origin, a 
behavior not found in the corresponding CTRW. We also 
show that certain CTRW solutions, converge extremely 
slowly toward the fractional diffusion approximation. 

The exact solution of the CTRW is based on a particu- 
lar choice of waitin g tim e and jump length distributions. 
Beginning in Sec. HI we investigate the fractional ap- 
proximation using a more general approach. Both cou- 
pled and decoupled CTRWs are considered. As far as I 
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know, the relation between coupled CTRWs and the frac- 
tional diffusion equation was not discussed previously. In 
Sec. ( V A ) we discuss Castiglione et al's objection to 
the fractional diffusion equation, for systems exhibiting 
strong anomalous diffusion. 



II. CTRW AN EXACT SOLUTION 

In this section we find an exact solution of the de- 
coupled, one dimensional, CTRW model in terms of an 
infinite sum of stable functions. Usually the solution of 
the CTRW, for finite times, is found using a numerical 
approach. The exact solution is used to understand the 
meaning and limitations of the fractional diffusion equa- 
tion. 

For the well known decoupled CTRW model, a particle 
is trapped on the origin for time ti, it is then displaced 
to xi = 5xi then the particle is trapped for time ^2 and 
then it jumps again, the process is then renewed. Let 
■4){t) be the probability density function (PDF) of the 
independent identically distributed (IID) random vari- 
ables {ti} i = 1, 2, • • •, while the IID displacements {5xi} 
are described by a PDF f[6x). The displacement 5xi 
is related to the coordinate of the particle according to 
5xi = Xi —Xi^i and xq = 0. Here it is assumed that start 
of observation is also start of the process. 

We assume that the PDF of waiting times -(/'(i) = 
la,i{t) is a one sided stable probability density. Namely, 
its Laplace t ^ u transform is ^(w) — exp(— u"), and 
< a < 1. The PDF of jump lengths is chosen to 
be a symmetric stable density f{5x) = 1^,q{2^/ ^6x), 
namely its Fourier transform is f{k) = exp(— |A;|^/2) and 
< < 2. The case fi — 2 (i.e., Gaussian jumps), 
< a < 1 and dimension d > 1 was considered in [10|. 
Properties of the stable densities can be found in [^8| - 31 1 . 
In the following sections the more general case where 
f{Sx) and ip{t) belong to the domain of attraction of the 
Levy stable laws, as well as coupled space time memories, 
is considered. 

Because the model is decoupled |l8] 



Pix,t) = Y,NcTis,t)Wix,s), 



s=0 



where NcTis,t) is the probability that s steps are made 
in time interval (0, t) and W{x, s) is the PDF that the 
random walk is on x after s steps. Using the convolution 
property of the symmetric stable densities, it is easy to 
show that 



Wix,s) 



(5) 



NcT{s,t) is found using the convolution theorem of 
Laplace transform (see details in [lO| ) 

NcT{0,t) = l-L„,iW, 



Nct{s, t) = Lc 



,1/q 



La,l 



(s + 1)1/"^ 



(6) 



and La,i{t) = la.i{t)dt is the one sided stable distri- 
bution function. Hence from Eqs. (§^) 

P{x,t) = [1- L^^i{t)]S{x)+ 



E 



La. 



t 



2I/M 



La 



M,0 



(s + l)i/c 
2^/^'x^ 



(7) 



The first term on the right hand side describes random 
walks for which the particle did not leave the origin 
within the observation time t; the other terms describe 
random walks where the number of steps is s. In Fig. 
(|l|) we show an exact solution of the CTRW process, in 
a scaling form, for the case a = 1/2 and /i = 1. 

The fractional diffusion approximation is reached when 
the summation in Eq. (^ is replaced with integration. 
We show below that such a replacement is not always 
valid. Using 



LaA { 1 /„ I —La.i 



1 t 

a s 



(s + l)i/a J ds VsV" 
t 



^ [ gl/a j ' 



(8) 




and neglecting the delta function contribution in Eq. (|^) 
we find 



P{x,t) 



(9) 

this approximation might be expected to work well only 
in the long time limit. The right hand side of Eq. (^) is 
the integral solution of the fractional diffusion equation 
(|l|) obtained by Saichev and Zaslavsky [i.e., only for 
< a < 1 and < /i < 2 and see Eq. (p4|) in Appendix 



(4) A] . In subsection H A we will show that in the vicinity of 



the origin x = and for fu, < 1 the fractional approxima- 
tion Eq. (^ does not work well. Let us therefore analyze 
the sum Eq. (|^) more carefully. 
First we rewrite Eq. (0) 



P{x,t)=y2W{x,s) / la.l{t')dt' 



(10) 



where W{x,0) = S{x). Using the Euler-Mclaurin sum- 
mation formula 1321 we have 



P(X, t) - d{x) I laAt')dt' + l-Wix, 1) 



la,l{t')dt' + 



2 



W{x,s) 



L,i (f) dt' 



ds 



(11) 



where • • • are the higher order terms in the Euler- 
Mclaurin formula. As discussed below, important cor- 
rection terms to the fractional diffusion approximation 
can be calculated based on Eq. ([ll]). In the long time 
limit the first two terms on the right hand side of Eq. 
( |rT| ) decay like t~". Within the fractional diffusion ap- 
proximation these terms are neglected. The third term 
on the right hand side of Eq. (|l|) yields the leading 
contribution to P{x, t) in the long time limit. For this 
term only contributions from large s are important, when 
t — > c». Using 



t/(s+l) 



we find in the limit < — > oo 



P{x,t) 



W{x,s)- 



as 



l/c 



,1/0 



ds. (13) 



This expression differs from the exact solution of the frac- 
tional diffusion equation by its nonzero lower limit in the 
integral. Comparing Eqs. and (|l^) we see that the 
shortcoming of the fractional approximation is that it 
attempts to give statistical weight to trajectories where 
number of jumps is "less than one" . Subtracting Eq. ( p^ ) 
from Eq. (^, and using la.iif) oc t~*^"+^) when i — > oo 
we have 



W{x,s)- 



t 



1/q+I 



as 



t 



,1/q 



ds (X t 



W{x, s)ds. 



The integral W{x, s)ds may become very large when 
x is small, and when x = the integral may diverge 
(i.e., since W{x, s)\x=o oc s~^/'^) . Then the convergence 
of the CTRW to the fractional approximation becomes 
extremely slow and when x = the fractional approxi- 
mation breaks down for fi < I. 




FIG. 1. We show t^^^P{x,t) versus the scaling variable 
^ = x'^/t for the CTRW process. The curves in tlie figure 
are for the times t — 1,5,25,125 (dots, dot dash, dashed, 
stared curve) respectively. We see that the CTRW solution 
converges in the limit of large times to the fractional diffusion 
approximation (solid line). The fractional solution diverges 
on the origin, hence it is cutoff in the figure. To find the 
exact solution we used Cauchy = 1 and Smirnov a = 1/2 
stable laws. The results shown in Figs. 1 — 4 were obtained 
using Mathematica. 



A. Diverging Solution of the Fractional Equation 

We now investigate in detail the behavior of the CTRW 
and the corresponding fractional diffusion equation at the 
origin x = 0. In Appendix A, the solution of the frac- 
tional diffusion equation is used to show that 



Pfrix,t)\,j:=o = 



2^*""/^ffi^'M.o(a;)U=o M > 1 



oo 



M< 1, 

(14) 

where lp,.o{x)\x=o = r(l//^)/(/^'''')- The subscript fr 
stands for fractional diffusion approximation. We as- 
sumed a < 1 since the case a = I yields stable propa- 
gator Pfrix,t) which does not diverge on the origin (see 
Appendix A). 
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In Appendix B the exact CTRW solution, Eq. (0), is 
used to find the behavior of the CTRW on the origin (i.e. 
for the non singular terms). For a < 1 we find 



^ ' ^' " I 7rr(l-Q) t" h<- — ^ 

(15) 

where C(z) is the Riemann zeta function and tpia) is the 
psi function [of course not related to ipit)]- In Fig- (H) 
we show P{x,t)\x=a versus t for fi — 1, exhibiting how 
the exact CTRW solution converges to its asymptotic 
limit at the origin. Comparing Eq. (15) with Eq. ( p^ ) 
we see that the infinity found for /i < 1 within the frac- 
tional framework is not related to the underlying CTRW. 
As mentioned, this shortcoming within the fractional ap- 
proximation is due to the fact that number of steps in the 
random walk is an integer which cannot generally be ap- 
proximated with a continuum approach [i.e., replacement 
of summation with integration in Eq. (^ is not justified 
for < /i < 1 in the vicinity of the origin]. In Fig. || we 
also show the approximation based on the Euler-Mclaurin 
formula Eq. ( |ll| ). In contrast to the fractional diffusion 
approximation, Eq. ( [ll| ) yields good agreement with the 
exact results. Note that divergence of the solution of 
the d > 1 dimensional fractional diffusion equation with 
/i = 2, a < 1 (Le. sub diffusive case) at the origin was 
discussed in |l0| , |3^ . For these cases the exact CTRW 
solution is a valuable tool. 



FIG. 2. The behavior of exact CTRW on the origin, 
P{x,t)\x=o versus t on log log plot for the case a = 1/2 and 
fi — 1. Not shown is the delta function contribution. The 
stars are the exact solution (^ while the solid curve is the 
asymptotic behavior Eq. (|l5|). The fractional diffusion equa- 
tion yields for this case P{x,t)\x^o = oo and hence is invalid. 
The approximation based on the Euler-Mclaurin formula, Eq. 
(prl), is the dot dashed curve. It yields good agreement with 
the exact results. 



B. Slow Convergence Toward Fractional 
Approximation 

Let us consider as an example the case a = 1/3 and 
jjL = 2. In Fig. ^ we show the exact CTRW solution 
in scaling form. As expected, for long times the CTRW 
solution seems to converge toward the solution of the 
fractional diffusion equation, though clear deviations of 
the CTRW solution from the fractional approximation 
are seen in the vicinity of the origin. In Fig. ^, a closer 
look at the behavior at the origin is presented. The figure 
shows that the CTRW convergence toward the fractional 
diffusion approximation is extremely slow, for t = 10^^ 
deviations from asymptotic behavior are still observed 
(note that since ip{u) = exp(— u^/'^) the natural time 
unit is 1, though the mean time between jumps diverges). 
Our improved approximation, Eq. (|ll|), yields a good de- 
scription of the underlying CTRW for intermediate and 
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long times. We note that as a ^ convergence of the 
CTRW solution toward the fractional approximation is 
expected to become much slower. And of course when 
a — !■ 1, convergence becomes faster, though then devia- 
tions from Gaussian behavior (i.e., /j, — 2, a = 1) become 
small. 



FIG. 4. The very slow convergence of the non singular 
CTRW solution toward the fractional approximation on the 
origin. We use the same parameters as in Fig. ^. The ex- 
act CTRW solution are the stars, the dot dash line is the 
fractional approximation, and the solid curve is the approxi- 
mation based on the Euler-Mclaurin formula, Eq. 




FIG. 3. The long time convergence of the exact CTRW 
solution toward the fractional diffusion approximation for 
a = 1/3 and ^j, — 2. We present t^^^ P{x,t) versus the scal- 
ing variable x^/t^^^. The CTRW curves in the figure are 
for the times t = 125,625,3125 (dot dash, dashed, dot) re- 
spectively, the solution of the fractional diffusion equation is 
the solid line. To obtain the exact CTRW solution we use 
the waiting time density /i/3,i(t) = 1/ {3TT)t~^^'^ K1/3 (^ffj) , 
where A'l/a is the modified Bessel function of the sec- 
ond kind. The jump length probability density is Gaus- 
sian f{5x) = (47r)~^/^exp(-5xV4). For a = 1/3 and 
fi — 2 the solution of the fractional diffusion equation is 
Pix,t) = 3|x|-^C/i/6,i(0 with C = t/\xf. 



III. THE GENERAL APPROACH 

While the exact solution presented in previous section 
gives insight into the validity of the fractional diffusion 
equation it is based on a particular choice of ijj{t) and 
f{Sx). Here we shall consider a more general approach. 

Let P{x,t) describe a normalized Green function of 
an unspecified one dimensional random walk; later we 
treat specific examples in some detail. We shall use the 
convention that the arguments in the parenthesis define 
the space we are working in, thus P(fc, u) is the Fourier- 
Laplace transform of the Green function P{x,t). It is 
assumed that P{k,u) is known exactly as is the case 
for different types of CTRWs and for generalized mas- 
ter equations of the type investigated in [|3j . 

Consider the expansion 



P(fc,u) = 5^A/„(m) 



{iky 



(16) 



where it is assumed that all the moments of the random 
walk 



M„{u) 



idk 



P{k,u) 



k=0 



(17) 



are finite. The case when the moments diverge will be dis- 
cussed later. According to Tauberian theorems Jl8| , p9t , 
the small u behavior of the moments in the Laplace do- 
main yield the long time behavior of the moments in the 
time domain. Namely, if A/„(u) ^ c„u"^"~^ when m — > 
then Mn{t) ^ c„t^"/r(l C„) for i ^ 00, where C„ > 
and c„ are constants. Our goal is to find the asymptotic 
moment generating function Pamg{^i'^) which contains 
all the information on the it ^ behavior of the mo- 
ments. We define this function according to 



PAMG{k,u) = I 



(18) 



ri=0 



If this series can be summed (see examples to follow) 
the function PAAicik^u) yields in a compact way all the 
asymptotic information on the moments of the under- 
lying transport process. If the inverse Laplace- Fourier 
transform of PAMG{k,u) is a normalized non-negative 
probability density, then it is safe to say that PAMcif^: 
is the Fourier-Laplace transform of the asymptotic Green 
function i.e., P{x,t) in the long time limit. This is the 
case for most Gaussian transport systems. Below we dis- 
cuss Levy walks where PAMG{k,u) does not yield the 
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Green function P{x, t) in the long time limit, even though 
it does contains all the information on the long time be- 
havior of the moments. 



IV. SUB-DIFFUSION 

We now consider as an example the well known decou- 
pled continuous time random walk in the sub diffusive 
regime. Similar to previous work (e.g. |l8|) we assume 



V'(w) = 1 - Au" + ciu^ 



< a < 1, 



(19) 



for u ^ 0, so that V'lO °^ for t oo; hence i\;{t) 

is moment-less. In what follows we use A — 1 restoring 
A only when it is important. We also assume 



f{k) 



2 



4! 



6! 



(20) 



where rm are the finite moments of the PDF f{5x). We 
assumed the f{Sx) is symmetric hence m2i+i — 0, the 
second moment being TO2 = 1. As discussed below some 
of these assumptions can be relaxed. 

The Green function P(x, t) of finding the random 
walker at x at time t is given in Fourier-Laplace space 
according to |l8| 



P{k, u) 



1 - '4>{u) 



1 



u l-i,{u)f{k)- 



(21) 



The long time behavior of this equation is usually investi- 
gated based on the long wave length approximation (e.g., 
namely by inserting 



Hu)f{k) = 1 - 



fc2 



(22) 



and ip{u) = 1 — in Eq. (pi|). This approach implicitly 
assumes that simple scaling Eq. (^) with /i = 2 holds. 
Let us now see why this is the case. 
We expand P{k,u) in k 

p(k, u) = - <^ 1 - n(u)— + \6n^(u) + m4n(u)] — — 

u [ 2 4! 



fc6 



(23) 



where f2(u) = ip{u)/[l — ijj^u)]. The i th term in the ex- 
pansion gives the moment Mi(u) of the random walker 
in-terms of the "microscopic" moments rrii and n{u). For 
example Afo = l/u which means that the normalization 
is conserved, M2{u) = VL{u)/u etc. In the long wave 
length approximation (or in the fractional diffusion equa- 
tion approach) one sets rrii = for z > 4. To see why 
and when this works well we must consider the high order 
moments Mi and i> A. For example 



Mi(u) 



[QVL^{u) + min{u)\ 



(24) 



We now consider the u — > limit of this expression since 
this limit will yield the asymptotic expression for M4^{t) 
when t ^ oo. Using Q.{u) ^ 1/u" it is easy to see that 



Miiu) 



6 



(25) 



and we see that Mji{u) is independent of 7714 in the limit 
u — > 0. Similar behavior is found for all the higher order 
moments 



M,,(.)^iiML, 



(26) 



which is independent of ci and m2j for j > \. This 
behavior is similar in some sense to normal (i.e., a = 1) 
random walk where all the moments Mn{u) converge in 
a limit to simple Gaussian behavior which is independent 
of the details of the underlying random walks. However 
now we are not considering Gaussian diffusion. It is also 
easy to show that behavior in Eq. (|^ is compatible 
with the scaling assumption Eq. (|). Using Eq. ( p^ ) one 
finds 



PAMG{k,u) = ^ 



fc^J(-l)J 

u(2u")J 



summing this geometric series we have 



P, 



AMG 



{k,u) 



,Q-1 



Or we may inverse Laplace transform Eq. (^ 
term and find 



PAAiG{k,t) = Y.[-T- 



r{l + aj) 



En 



(27) 

(28) 
term by 



(29) 

where Ea{x) is the Mittag-LefHer function. Since Eq. 
( p8| ) is the Fourier Laplace transform of a non nega- 
tive probability density (see e.g. Appendix A), it yields 
P{x, t) in the long time limit. 

Eq. (|2^) is the Fourier-Laplace transform of the frac- 
tional diffusion equation (|^) when /i = 2 and < a < 1, 
and well known within the CTRW community |^ . The 
inverse Fourier-Laplace transform of Eq. ( p8| ) was inves- 
tigated in HJl^,^ in dimensions d — 1,2,3. Here we 
showed that: {i) this equation does indeed describe the 
long time behavior of the moments of the random walk 
to all orders and (ii) that these moments depend only 
on three parameters of the model a, m2 = 1 and A ~ 1 
(i.e., universality) . Our approach clarifies the usual long 
wave length approximation which is based on the exact 
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calculation of only the first two moments. In the follow- 
ing section we will discuss coupled memories where the 
asymptotic behavior is not as straightforward as for the 
decoupled case. 

In our derivation we assumed that start of observation 
and start of the process coincide. If the first step is de- 
scribed by tpi{t) ^ ip{t), one can show that our results are 
still valid for i{ji{t) decaying faster then tpi(t) oc t~^^^^\ 
with P > 0. When (3 << a the convergence becomes 
slow. If the random walk is biased, mi ^ 0, one can 
easily show that biased fractional diffusion equation (e.g. 
[pO|) holds in the long time limit. 

V. ENHANCED-DIFFUSION 

We now consider an example exhibiting enhanced, 
Levy walk type of diffusion. We start by introducing the 
coupled CTRW jump model, investigated by Zumofen, 
Klafter and Shlesinger [^,0 in the context of chaotic 
maps. Such a random walk is also related to transport in 
random media ]3^ , |39| ], tracer diffusion in turbulent flow 
pof and to the blinking of Quantum dots [^. Closely 
related models are the velocity models investigated in 

In CTRW the random walk is entirely specified by 
ip{Sx,t), the probability density to move a distance Sx 
in time i in a single jump event. For the jump model a 
coupled space-time memory is assumed 



4>{Sx,t)^^S{\5x\-t)ij{t)- 



(30) 



Such a model describes a particle trapped on the origin 
for time ti then it jumps to a new location whose distance 
from the origin is |fe| = ti (i.e, xi = zLti with equal prob- 
abihty), then the process is renewed. From Eq. ( |30| ) we 
see that the random times {ti} are distributed according 
to tl}{t) and the length of each jump i is \Sxi\ — ti. Hence 
a large jump will "cost" a long time. This is different 
from the Levy flight model where jumps on all scales are 
performed at constant time intervals. Thus in some ap- 
plications Levy walks are considered more physical than 
Levy flights; however, as discussed below these two mod- 
els are in fact deeply related. The space time coupling 
in Eq. (|3^) guarantees that for the Levy walk model, 
P{x, t) = for \x\ > t, this in turn implies that all mo- 
ments of the random walk are finite. For Levy flights 
even moments diverge. 

The Laplace-Fourier transform of the Levy walk Green 
function is [56|l 



P{k,u) = 



1 - V(u) 



u 1 — ■0 (u, fc) 



ip{k,u)^ / dte^"* cos {kt)ip{t). 
Jo 



(32) 



The moments of the random walk are now calculated 
using Eq. ( ^ and Mathematica, one finds Mo{u) = l/u 



M2{u) 



1 7/^(2) (y) 



M^u) = - 



u 1 — tp{u) 
1 [l-i^(u)]i/;(4)(u) + 6[V'(2) (u)]'- 



(33) 



where ip^^^\u) is the 2j th derivative of ip{u) with re- 
spect to u. Odd moments vanish due to the assumed 
symmetry of the random walk. Higher order moments 
are calculated in a similar way, for the sake of space they 
are not included here. 



A. Sub— Ballistic Enhanced Diffusion 

Let as now consider ip{u) ^ 1 ~ ut + Bu^ ■ ■ ■ and 
1 < (3 < 2. Unlike the previous example, now the mean 
waiting time r is finite and the second moment of the 
waiting time distribution diverges. The model exhibits 
enhanced diffusion (x^) cx and a — i — (3. 

As we shall see in detail, the model exhibits a strong 
type of anomalous diffusion. By definition |27| strong 



anomalous diffusion behavior exhibits M2j{t) cx t^^^^ 
where /(? ') is a non-linear function of j (see related work 
48 1). Castiglionc et al |27| point out that strong 
anomalous diffusion implies the failure of the standard 
scaling assumption Eq. (^), since this equation predicts 
M2j {t) (x t"^ a behavior called weak anomalous diffusion. 

Castiglione et al argue quite generally that a dynami- 
cal system exhibiting strong anomalous diffusion cannot 
be described by fractional diffusion equation. We show 
below, based on work of Zumofen et al ||2j] and others, 
that the Green function P{x, t) is well described by the 
fractional diffusion approximation. 

First we consider a long wave length approximation, 
and show the relation between this approximation and 
the fractional diffusion framework. We rewrite Eq. ( |3l| ) 
in the form 



uP{k, u) — wp{k, u)P{k, u) = 1 — ip{u). 



(34) 



Due to Tauberian theorems, the behavior of P{x, t) for 
t ^ CX) is controlled by the behavior of P(fc, u) at u — > 0, 
hence in the small u limit (and fixed k) one finds 



(31) P(fc,u)il 



where 



%l){k,u)\u=o — u j t cos (kt) ip{t)dt + ■ 
Jo 

T - Bu^-^ + ■■■. (35) 
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We now consider the fc ^ limit using il;{k,u)\u=o 
1 + B\k\l^ cos (7r/3/2) + • • • leading to 



P{k,u) 



1 



u + B\k\P\cos {-kPI^) \/t 



(36) 



Eq. ( |3q ) is rewritten in terms of a fractional diffusion 
equation using convenient units 



dPfrix,t) 1 



dt 



2d\x\ 



■Pfr{x,t) = 5{x)5{t). 



(37) 



Eq. ( pT) ) describes Levy flights, whose solution is a sym- 
metric Levy stable PDF given in Appendix A, Eq. (5f ). 
The derivation of Eq. (^^ based on the long wave length 
approximation is not rigorous; however Zumofen et al 
used a numerical inverse Fourier-Laplace technique 
to show that the propagator is well described by a Levy 
stable PDF. This result was verified by several authors, 
Araujo et al [^3| used a numerically exact enumeration 
technique ji^ and Mantegna and Weron and Weron 
p6| used a Monte Carlo approach. As briefly mentioned 
in the introduction Kotulsky used a rigorous limit 
theorem approach to reach the same conclusion. Thus 
in contradiction to the claim made in [ pTj , the fractional 
equation yields a meaningful approximation to the un- 
derlying strong anomalous diffusion process under inves- 
tigation. 

We are still left with a puzzle: the fractional equa- 
tion predicts a non analytical behavior oi P(k, u), namely 
the divergence of the even moments of the random walk, 
while we know that these moments, for any finite time 
are finite. Let us therefore investigate the moments in 
greater detail, using Eq. (^) Mo(u) — 1/u and for j > 1, 



T 



/3-2-2J 



(38) 



where 62j(/3) = ^iLq^\^ ~ P\- Inverting to the time do- 
main we find the mentioned strong type of anomalous 
diffusion M2j{t) cx t^i+i-/? for j = 1, 2 • • • while Mq = 1. 
We now investigate the behavior of the asymptotic mo- 
ment generating function, based on the method in Sec. 
n|. According to Eq. (|l|) 



using the identity (obtained using Mathematica) 



, (39) 



= 1 



we find 



(2j)! 



-1 + (l + x)'^/'^ cos [/3arctan (Vx)] 

(40) 



PAMG{k,u) 



B 



,/3-2 



(41) 



We note that unlike the sub diffusive case, PAMG{k,u) 
Eq. ( pl| ) is not the Fourier Laplace transform of the 
asymptotic P{x, t) since limj.^00 PAMci^^ ") = when 
u is fixed. 

To conclude, according to long-wave length approxi- 
mation and previous work, the jump model propagator 
P{x, t) is described by the fractional diffusion equation 
with (a — = j3). This approximation does not de- 
scribe the behavior of the moments of the underlying 
random walk including the second. These are described 
by the moment generating function Eq. (^l|); thus, two 
functions yield the details on the long time behavior of 
the underlying random walk. And strong anomalous dif- 
fusion docs not necessarily imply the breakdown of the 
fractional approximation, though one should take care in 
the interpretation of the results obtained by it. 

A similar situation occurs in the field of inhomogeneous 
line broadening [ pO[ , where the line is well approximated 
by Levy stable densities |5^] (due to long range interac- 
tions between defects and chromophores). However (due 
to cutoffs) even moments of the line exist [Q. As dis- 
cussed by Stoneham these moments are sensitive to 
behavior of the line in its wings, and hence in the usual 
experimental situation (in the field of line broadening) 
are not considered relevant. 



B. Ballistic Diffusion 

We now briefly consider the jump model with ip{u) = 
1 — u'^ • • •, with < /? < 1. For this case the model ex- 
hibits ballistic diffusion {x'^(t)) cx [Q. Without going 
into details, we find the asymptotic moment generating 
function 



PAMcik^u) 



[u- 



ikf + {u-tkf' 



(42) 



Zumofen et al P4| used an expansion in the small pa- 
rameter u ± ik and numerical simulation for f3 = 1/2, 
showing that Eq. (|4^) describes well long time behav- 
ior of P{x,t) (at least for (3 — 1/2). Here we showed 
that the approximation in [p4| , yields the behavior of the 
moments to infinite order. Eq. ( ^2| ) is not related to a 
known fractional diffusion equation. Note that the frac- 
tional diffusion equation in the ballistic limit = 2, a = 2 
yields the wave equation. 



VI. LEVY FLIGHTS WITH LONG RESTS 

We now consider decoupled CTRW, where ip{Sx, t) = 
^j{t)f{5x), assuming the moments of f{Sx) diverge. 
Clearly, now the general approach presented in Sec. Ill 
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breaks down and a generalization is now considered. We 
assume 



fik) = 1 - 



4! 



4! 



(43) 



and < /i < 2. The coefficients called the 

amplitudes of the PDF f{Sx). An example being the 
symmetric stable densities f{k) = Z^_o(fc/\/2) where 
lfj.fi{k) = exp(— |fc|'^). We assume as before that 



1p{u) = 1 - + CiM^" 

for u — > 0, then using Eq. (El 



< a < 1, 



(44) 



1 f \k\^^ IfcP^ 

p{k, u)^-li- n{uy-^ + [6n\u) + a2^n{u)] 



\k\^ 



(45) 



This equation has a structure similar to Eq. (|2^), we 
see that the amplitudes aj^ are natural generalizations 
of the moments m2j- We define the amplitudes Ajf^{t) 
according to 



P{k,u) 



3=0 



it follows from the discussion in Sec. IV that in the limit 
u ^ 0, Aj^{u) - (2j)!/(2%J°+i). We define the am- 



plitude generating function, in the spirit of Eq. (18) 
according to 



PAG{k,u) = Y,i-^y 

and it is easy to show that 



|fc| 



2JuJ'= 



(47) 



(48) 



The right hand side of Eq. (48) is the Fourier-Laplace 
transform of the fractional diffusion equation (|l]) . Hence 
the fractional equation describes the amplitudes of the 
random walk in the limit i — > oo. However this does not 
necessarily imply that for all x the corresponding P{x, t) 
describes the CTRW in the limit i — > cx). As shown al- 
ready in subsection HA, for /i < 1 and on the origin 



where P{x, t) attains its maximum, the CTRW solution 
and the fractional diffusion equation solution arc differ- 
ent. This limitation of the fractional equation is related 
to the fact that it is based on a small k expansion which 
implies large x behavior. Note that Eq. (^ ) was recently 
suggested by Kutner |52| to describe a Weierstrass flight. 



VII. DISCUSSION 

The fractional diffusion equation yields the asymptotic 
behavior of coupled and decoupled CTRWs. However 
the fractional approach has its limitations if compared 
with ordinary diffusion approximation. Careful analy- 
sis of the underlying random walk is needed for a bet- 
ter understanding of the domain of validity of the frac- 
tional equation. In particular the fractional approxima- 
tion (for decoupled CTRWs) breaks down at the origin 
X = (for {0 < a < 1} n {0 < /I < 1}). More gener- 
ally the convergence of the CTRW solution, at a: = 0, 
toward the fractional approximation may become ex- 
tremely slow. This behavior is different from ordinary 
random walks where (usually) (i) deviation from Gaus- 
sian behavior is found in the tails |x| oo (i.e., for long 
though finite times) and (ii) convergence toward the frac- 
tional approximation is typically fast. Also note that 
the fractional equation has a CTRW foundation only 
when {0 < a < 1} n {0 < /Lt < 2} while the regime 
{1 < a < ^1 < 2} is not related to the CTRWs under in- 
vestigation (see however work in ^ for {1 < a < 2 = /x}). 

If they exist, moments of the decoupled CTRW con- 
verge in a long time limit to the behavior predicted by 
the fractional equation. For coupled memories describ- 
ing Levy walks the situation is more complicated. The 
fractional approximation describes the asymptotic long 
time behavior of the Green function P{x,t), though it 
does not describe correctly the moments of the underly- 
ing random walk, not even the second moment. In this 
case we may characterize the random walk using both 
P(x, t) and the asymptotic moment generating function. 
These two function yield different types of information, 
it is still to be seen if the asymptotic moment generating 
function has any universal features. 

We note that a similar situation exists also for some 
simple random walks which are approximated with the 
ordinary diffusion equation. To see this consider as an 
example the sum of A'^ oo independent, identically 
distributed random variables {xi}, i = 1,2---A^. As 



well known the sum ^ 



N 



converges in a limit to a 



Gaussian behavior, provided the variance of Xi exists. 
Assume the variance exists but higher order moments of 
Xi diverge. Then clearly the central limit theorem, or 
identically the ordinary diffusion approximation, fails to 
predict correctly the behavior of the high order moments 
of the random walk (i.e., the Gaussian central limit the- 
orem does not hold at the tails of the Green function). 
The situation for the fractional diffusion equation is sim- 
ilar to this case, in that it fails to predict correctly the 
behavior of the moments of the Levy walk (i.e., the Levy 
central limit theorem does not hold at the tails of the 
Green function of the Levy walk, where P{x, t) — for 
\x\>t). 

Note that our conclusions are valid only for free bound- 
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ary conditions, for other boundary conditions we know 
little on domain of validity of the fractional diffusion 
equation especially when fi < 2 (see however work in 

Sill)- 
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(54) 



this equation being valid for < a < 1 and < /i < 2. 



IX. APPENDIX B 



VIII. APPENDIX A 

We investigate the solution of fractional diffusion equa- 
tion, with a < 1 and < jj, < 2. Following ||] we rewrite 
the solution in Fourier-Laplace space 

/"OC 

Pfr{k,u)^u'-'^ dse-^("°+l'=l"/2). (49) 



Since the symmetric Levy stable probability density 
lp,,a{x) and exp(— |fc|'^) are Fourier pairs 

Pfr{x,u) = 2^'t^u^-^ dse-'^^l^fl \-s^) 

(50) 

For a = 1 we use the Laplace pair exp(— su) and S{t — s) 
and find as expected 



We investigate the CTRW solution 



1 



u 



-Six) 



5 = 1 



(55) 



The first term can be easily handled yielding La^i{t)S{x) 
decaying for long times like S{x)t~°'. This singular term 
is neglected in the fractional diffusion approximation. 
Omitting this term we find for x = 0, 



P{x,u%^o 



1 - 



1/m 



Ifj^M (x) \x=a 



(51) 



and when fi = 2 the solution is Gaussian. 

We now consider a < 1 and investigate the behavior 
on the origin. Using Eq. (|50| ) 



(52) 



where Z^,o(a;)U=o = r(l/Ai)/(Ai7i"), we find 



2i/p ( 1 — I poiyLog 



-, e 



lp.,0 (x) U=o (56) 



where PoiyLog [n, z] is the n th polylogarithm function 
of z. For fi < 1 and u — s- we find 

P{x, w)U=o - 2i/^u"-iC lM^)U-o, (57) 

where (^(z) is Riemann's zeta function. For /i = 1 we use 
^exp(— sw")/s = — ln[l — exp(— su")] to find 



Pfr{x,u)\x^o 



2i/MWp-ir (1 _ 1/^) (x) > 1 



oo 



(53) 



1 



P(x,u)U=o - -2u"-Mn(u")Zi,o(x)U=o- (58) 

For 1 < /i < 2 we use the Euler-Mclaurin summation 
formula 



which when inverted yields Eq. ([l^). For < 1, 
Pfrix,u)\x=o is infinite due to the zero lower bound in 
the integral Eq. (|2|) . 

An integral solution of the fractional diffusion equa- 
tion is found using the Laplace pair exp(— su") and 
{l/a){t/s^+''^°'%,oit/s^/°') where /„,i(i) is the one sided 
stable density whose Laplace transform is exp(— m"). 
Hence 

Pfrix,t) = 



-su" -1//J 



dk 



l-|-^e-'^"°(/c + l)-i 



fe=i 



1 Id er^" 



+ 



(l + fc)i/P 2 12 dfc (1 + A:)^/^ ° 

(59) 



where • • • are the higher order terms in the Euler- 
Mclaurin formula. The integral in Eq. ( p9[ ) in the limit 
M — > yields 
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poo —ku^ 

/ dk— — — 

lo (1 + A:)^/'^ 



/"/'^-"r - ^) ' (60) 



provided that /i > 1. This term is much larger than the 
other terms in Eq. (59); hence, we find 



p{x, u)i^o - 2i/^W^-ir (^1 - lMU=o- (61) 

Inverting Eqs. (I^HJe^) we find Eq. (|l|). 
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